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The  kinetics  of  structured  relaxation  in  a  two-dimensional  model  atomic  glass 
quenched  infinitely  rapidly  from  the  melt  to  0.55  of  the  glass  transition  tempera¬ 
ture  was  simulated  by  the  molecular  dynamics  method  to  study  the  chronological 
ordering  of  the  atomic  kinematics  associated  with  such  relaxations.  Over  the 
very  short  periods  of  aging  (c.a.  200  atomic  fluctuations)  accessible  to  the  MD 
method,  a  Williams- Watts  form  of  relaxation  with  a  fractional  exponent  of  0.5 
was  found  to  hold  for  excess  enthalpy,  free  volume,  and  site  distortion  parame¬ 
ter.  The  distribution  of  free  energy  barriers  associated  with  this  relaxation  that 

resulted  from  the  analysis  could  be  scaled  up  to  describe  processes  occuring  on 

r  * 

macroscopic  time  scales,  and  agrees  well  with  experimental  results  in  CuzZrx_z 
glasses.  Results  on  the  clustering  of  relaxations  and  other  topological  features  of 

A 

the  relaxation  process  are  also  reported.  t<-  f  ^  M  Cc' 
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I.  INTRODUCTION 


The  atomic  details  of  structural  alterations  of  amorphous  solids  during  aging 


below  their  glass  transition  temperatures  and  the  consequences  of  this  on  sub¬ 


sequent  inelastic  behavior  has  been  of  considerable  interest.  Simulation  of  such 


complex  atomic  processes  by  computer  molecular  dynamics  (MD)  has  been  very 


informative  in  the  past.  Among  the  many  such  simulations,  those  of  Takeuchi 


and  coworkers  [1,2]  on  two  componer.c  atomic  glasses  of  the  CuIZr1_I  type, 


carried  out  in  three  dimensions  have  resulted  in  considerable  insight  on  struc¬ 


tural  relaxations.  In  spite  of  this,  however,  both  the  kinematics  and  the  ki¬ 


netics  of  such  relaxations  remain  inadequately  understood.  We  have  carried 


out  a  detailed  MD  simulation  using  the  same  truncated  CuZr  pair  potentials 


used  by  Kobayashi,  et  al  [3],  but  in  two  dimensions  for  the  specific  purpose  of  a 


more  thorough  understanding  of  the  atomic  mechanisms  of  plastic  flow  in  simple 


atomic  amorphous  media.  In  three  related  communications,  we  have  reported 


results  of  a  simulation:  on  the  melting  and  glass  transition  process  in  a  two¬ 


dimensional  idealized  material  [4]  (referred  to  here  as  (I));  on  the  topological 


features  of  structural  relaxations  in  such  a  medium,  both  in  subcooled  melts 


and  in  solids  below  their  Tg  [5]  (referred  to  here  at  II);  and  on  the  details  of 
large  strain  plastic  shear  [6]  (referred  to  here  as  IV).  Here,  in  this  same  series  of 


studies,  we  report  on  the  kinetics  of  structural  relaxations  below  Tg  in  the  same 


two-component  amorphous  medium  modeled  in  two  dimensions.  We  also  give 


additional  details  of  some  of  the  topological  features  in  the  material,  which  is 


more  prone  to  undergo  structural  relaxation  than  the  average. 
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II.  DETAILS  OF  THE  SIMULATION 


2.1  The  Simulation  Cell 

The  details  of  the  MD  simulation,  the  choice  of  the  interatomic  pair  poten¬ 
tial  of  interaction,  its  truncation,  the  two-  dimensional  simulation  cell  with  its 
periodic  boundary  conditions,  and  the  normalization  of  parameters  were  given 
in  (I).  Additional  details  on  the  intensive  state  parameters,  such  as  atomic  level 
stresses,  moduli,  enthalpy,  site  distortion,  and  local  free  volume  were  given  in 
(II),  together  with  many  results  on  the  topological  details  of  structural  relax¬ 
ations.  The  results  of  the  simulation  we  present  here  were  on  the  same  rect¬ 
angular  cell  containing  144  atoms  of  two  types  representing  Cu  and  Zr,  placed 
initially  randomly,  and  melted  and  quenched,  as  described  in  (I).  The  simula¬ 
tions  reported  here  were  at  a  temperature  of  T*  =  0.1  (=  0.4T*m)  and  under  a 
constant  external  normalized  pressure  of  p  =  1.0.  This  temperature  is  well  below 
the  glass  transition  temperature  of  T* i  =  0.18.  The  starting  state  of  the  mate¬ 
rial  was  one  which  had  seen  an  essentially  infinite  quenching  rate  (path  EF  in 
Fig.  1  of  II)  from  the  meliing  temperature  of  0.25,  i.e.,  the  starting  material  had 
a  Active  temperature  equal  to  the  equilibrium  melting  temperature.  Many  topo¬ 
logical  details  of  structural  relaxation  of  this  quenched  glass  were  reported  in  II. 
The  magnitude  of  the  so-called  boundary  mass  used  in  the  present  simulation  to 
maintain  the  external  pressure  constant  was  4.  The  concept  of  a  boundary  mass 
associated  with  the  volume  of  the  simulation  cell  as  an  extra  degree  of  freedom 
to  maintain  the  external  pressure  constant  was  introduced  by  Andersen  [7].  It 
has  been  widely  used  in  MD  simulations  where  only  equilibrium  structural  prop- 


erties  are  of  interest,  and  the  kinetics  of  reaching  the  equilibrium  is  not.  In  this 
simulation,  however,  where  the  kinetics  of  structural  relaxation  were  of  equal 
interest  as  the  topological  details  of  the  atom  exchanges,  the  proper  choice  of 
the  boundary  mass  was  an  important  consideration. 


2.2  Choice  of  Boundary  Mass  in  the  Simulation 


In  the  method  of  Andersen  [7],  simulation  under  constant  external  pressure 
is  accomplished  by  the  introduction  of  a  new  dynamic  variable  Q  that  represents 
the  changes  in  the  overall  volume  of  the  simulation  cell.  Then,  the  equations  of 
motion  of  all  n  atoms  in  the  system  and  the  cell  volume  are  derivable  from  a 
modified  Lagrangian  L  given  by  1 


L(Pi,Pi,Q,Q)  =  !W/3X>  •  Pi  -  l't4>[Ql/*Pij)  +  \mQ2  -  aQ,  (1) 


where  p,  =  (r./Q1/3)  are  dimensionless  position  coordinates  of  the  n  atoms,  p,, 
their  velocities,  p,; ,  the  dimensionless  relative  separations  of  atoms  i  and  j,  <f> 
the  pair  potential,  m  the  mass  of  an  atom,  Q  the  coordinate  representing  the 
volume  viewed  as  a  one  dimensional  entity,  such  as  the  length  of  a  cylinder 
capped  by  a  piston,  a  a  constant  to  be  interpreted  as  the  external  force  result¬ 
ing  from  the  external  pressure  acting  across  the  boundaries  of  the  simulation  i  Fur 

I 

cell  in  a  work  conjugate  sense  with  the  volume  coordinate  Q,  and  finally  M,  the 

1  All  symbols  not  specifically  defined  here  are  the  same  that  have  been  used  in  (I)  [4],  where 
a  lull  list  of  symbols  can  be  found. 
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so-called  boundary  mass.  In  Eqn.  (1),  the  first  two  terms  on  the  RHS  repre¬ 
sent  the  normal  Lagrangian  of  the  system  of  n  particles  if  Q  is  interpreted  as 
the  volume.  The  third  and  fourth  terms  relate  to  the  Lagrangian  of  the  total 
volume  of  the  cell  in  the  interpretation  of  Q  as  the  added  dynamic  variable  rep¬ 
resenting  the  volume.  The  generalized  velocities  and  momenta  that  result  from 
this  scaled  Lagrangian,  and  how  they  can  be  related  to  the  real  system  has  been 
discussed  in  detail  by  Andersen  [7]  and  need  not  be  repeated  here.  It  should 
suffice,  however,  to  observe  that  the  introduction  into  the  system  of  a  pseudo 
coordinate  to  represent  the  volume  and  a  pseudo  mass  associated  with  it  permits 
maintaining  the  external  pressure  constant  by  permitting  the  volume  10  change 
as  necessary.  Clearly,  however,  this  introduces  into  the  simulation  system  a  new 
normal  mode  and  a  new  momentum  with  an  added  normal  mode  frequency  that 
can  interact  with  the  other  physically  meaningful  modes  of  atoms,  and  may  re¬ 
sult  in  unwanted  resonances.  It  can  be  shown,  and  is  also  readily  clear  that  the 
introduction  of  this  pseudo  coordinate  does  not  affect  the  equilibrium  properties 
of  the  system  and  its  topological  features  in  equilibrium.  It  will,  however,  influ¬ 
ence  the  kinetics  of  reaching  equilibrium,  and  the  topological  features  of  non¬ 
equilibrium  states  that  are  of  interest  in  the  simulation.  Andersen  has  proposed 
that  unwanted  artifacts  in  the  kinetics  can  be  reduced  if  the  boundary  mass 
M  is  chosen,  so  that  the  period  of  the  fluctuation  of  Q  in  the  scaled  system  is 
approximately  equal  to  Q1/3  divided  by  the  speed  of  sound  in  the  simulation 
medium,  when  Q  is  equated  to  the  volume.  For  a  simulation  cell  containing 
144  atoms,  this  rule  results  in  a  boundary  mass  M  =  3.65 m  of  individual  atom 
masses.  We  have  taken  4  as  the  nearest  integer  value  to  be  used  in  the  simula- 
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tion.  In  addition,  however,  we  have  tried  simulations  of  structural  relaxations 
with  boundary  masses  of  400  and  4  x  10~2  to  determine  the  sensitivity  of  the 
results  on  this  choice.  Thus,  as  presented  in  Appendix  I,  while  a  choice  of  400  m 
for  the  boundary  mass  produced  overdamped  conditions,  very  markedly  slowing 
down  the  structural  relaxation,  a  choice  of  4  x  10-2  produced  erratic  results. 

III.  RESULTS 

3.1  Time  Dependent  Changes  in  Structural  Properties 

The  time  dependent  changes  in  several  ensemble  average  intensive  structural 
properties  for  a  period  of  200  atomic  fluctuations  (4000  time  steps)  are  shown  in 
Figs,  la-lg.  Here,  the  period  of  an  atomic  fluctuation  is  taken  as  5.4xl0lssec. 
(see  (I),  [4]).  The  figures  show  both  the  level  of  fluctuations  in  overall  average 
properties,  as  well  as  their  monotonic  changes  over  this  period,  first,  for  the  en¬ 
tire  ensemble  of  144  atoms,  but  also  for  the  subensembles  of  Zr  and  Cu  atoms, 
as  marked.  Among  these,  Fig.  la  shows  the  changes  in  the  average  Voronoi 
polygon  volume  per  atom  (r2  times  area  per  atom).  Clearly,  the  regular  fluctu¬ 
ations  due  to  the  boundary  mass  are  apparent  -  not  only  in  the  fluctuations  of 
the  volume,  but  to  some  extent  in  all  other  properties,  as  can  be  seen  by  inspec¬ 
tion  of  Figs,  lb-lg.  Figures  lb  and  lc,  e.g.,  show  the  corresponding  changes 
in  the  ensemble  average  values  of  the  enthalpy  and  site  distortion  parameter  2, 
while  Figs,  ld-lg  give  the  related  changes  in  the  atomic  site  pressure,  maximum 
shear  stress  (in  the  plane  regardless  of  orientation),  the  xy  shear  stress,  and  the 

aFor  definitions  of  these  intensive  parameters  describing  structural  properties  of  the  glassy 
state,  refer  to  (II)  (5 j. 


atomic  site  bulk  modulus  -  all  averaged  over  the  entire  ensemble  as  well  as  over 
only  the  Zr  atoms  and  the  Cu  atoms,  as  indicated  in  the  figures.  The  averaging 
has  been  done  for  every  fluctuation,  i.e.,  after  every  20  time  steps  without  any 
smoothing.  Clearly,  as  stated  above,  the  major  and  regular  fluctuations  with  a 
period  of  about  18  atomic  fluctuations  is  due  to  the  border  mass  and  must  be 
ignored.  When  these  are  discounted,  there  is  a  clear  trend  of  decrease  in  the 
volume,  enthalpy,  and  site  distortion.  Large  local  fluctuations,  however,  remain 
particularly  in  the  atomic  level  stresses  and  local  bulk  moduli.  Additionally,  it  is 
evident  that  the  atomic  site  pressures,  maximum  shear  stresses,  and  bulk  mod¬ 
uli  for  the  Zr  atoms,  are  much  larger  than  the  average,  while  the  reverse  is  true 
for  the  Cu  atoms.  As  already  noted  in  (II),  this  is  a  direct  result  of  the  stiffer 
environment  of  the  Zr  atom  than  that  of  the  Cu  atom.  The  effect  of  this  on  the 
enthalpy  fluctuations  is  in  the  reverse  order,  where  as  must  be  the  case,  more 
energy  is  stored,  on  the  average,  in  the  more  compliant  surroundings  of  the  Cu 
atoms.  Furthermore,  it  can  also  be  observed  that  in  the  regular  oscillations  of 
the  boundary  mass,  the  fluctuations  of  the  Voronoi  volume  and  the  atomic  pres¬ 
sure  are  out  of  phase,  i.e.,  as  the  atomic  site  pressure  increases,  the  local  volume 
decreases,  and  the  atomic  site  bulk  modulus  increases.  This  direct  response  of 
the  local  site  properties  under  the  action  of  the  boundary  mass,  however,  does 
affect  the  kinetics  of  structural  relaxation  by  speeding  it  up  to  some  extent  (see 
Appendix  I).  We  note  from  Figs,  la-lg  that  there  are  substantial  reductions 
in  excess  properties  (over  and  above  those  of  an  ordered  hexagonal  crystal)  in 
enthalpy,  local  volume  per  atom  (Voronoi  volume),  and  site  distortion  over  the 
total  relaxation  span  of  200  fluctuations.  There  is  a  corresponding  increase  in 
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the  local  bulk  moduli  from  17  to  about  18  (in  units  of  dimensionless  pressure) 
or  by  an  amount  of  6%.  In  comparison,  there  are  no  important  relaxations  in 
the  atomic  level  pressures  and  maximum  shear  stresses.  Since  the  changes  in 
enthalpy  must  come  from  the  reduction  in  root  mean  square  elastic  strain  en¬ 
ergies,  however,  it  must  be  concluded  that  the  observed  changes  come  about 
primarily  from  the  increases  in  bulk  modulus  and  a  corresponding  increase  in 
atomic  level  shear  moduli.  Although  all  ensemble  average  properties  fluctuate 
considerably,  the  ratio  of  the  average  fluctuation  amplitude  to  the  current  av¬ 
erage  level  of  the  property  shows  large  differences  between  different  properties. 
These  ratios  are  listed  in  Table  I  and  indicate  that  initially,  the  ratios  are  largest 
for  the  atomic  pressure  and  shear  stress,  while  they  are  comparatively  low  for 
the  atomic  volume.  Over  the  simulation  span,  however,  the  fractional  changes 
in  the  atomic  volume  exceed  those  in  the  stresses.  It  is  also  interesting  to  note 
that  although  the  change  in  the  average  levels  of  the  atomic  stresses  are  quite 
small,  they  are  very  substantial  in  the  amplitudes  of  fluctuation.  The  source 
of  these  differences  are  clear.  Thus,  while  the  fluctuations  and  changes  in  the 
atomic  volume  are  quite  small  on  the  absolute  scale,  being  limited  by  the  overall 
level  of  free  volume  that  can  be  squeezed  out,  the  stresses  can  fluctuate  between 
large  negative  to  large  positive  values  about  an  equilibrium  level  in  response  to 
only  small  displacement  amplitudes.  To  concentrate  attention  more  specifically 
on  the  actual  fluctuations,  unaffected  by  the  superimposed  large  fluctuations 
produced  by  the  boundary  masses,  the  root  mean  square  deviations  from  the 
current  ensemble  average  values  in  the  structural  parameters  were  determined 
after  each  atomic  fluctuation  (for  every  20  time  steps)  according  to  the  formula 
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xm,  =<  <  X(i)  >t)2  >)'2  . 


(2) 


In  Eqn.  (2),  the  brackets  mean  ensemble  averages  either  over  the  entire  pop¬ 
ulation  or  over  only  the  Zr  and  Cu  atoms,  and  t  indicates  that  the  information  is 
evaluated  at  a  given  time  (after  every  20  time  steps).  The  results  of  this  type  of 
evaluation  are  given  in  Figs.  2a-2g,  paralleling  the  results  given  in  Fig.  1.  These 
results  show  some  of  the  same  effects  shown  in  Fig.  1.  The  fluctuations  now  do 
not  show  the  overwhelming  pulsation  resulting  from  the  boundary  mass.  There 
is  a  clear  and  perceptible  decrease  in  the  fluctuations  in  the  atomic  volume  and 
the  site  distortion  parameter  over  the  simulation  span.  Since  the  average  levels 
for  the  fluctuations  of  the  Zr  atoms  differs  markedly  from  that  of  the  Cu  atoms, 
the  root  mean  square  fluctuations  in  either  of  these  atoms  is  less  than  that  over 
the  entire  ensemble  for  most  of  the  properties,  by  virtue  of  the  mode  of  the 
overall  ensemble  averaging  procedure. 


Much  experimental  research  on  structural  relaxation  in  amorphous  media, 
in  common  with  stress  relaxation  and  internal  friction,  has  established  that 
such  relaxations  are  not  simple  mono-energetic  processes  characterizable  by  a 
single  relaxation  time.  It  has  instead  been  proposed  by  many  researchers  that 
the  relaxations  in  an  amorphous  solid  are  distributed  in  relaxation  time,  or  in 
a  somewhat  more  fundamental  variant,  that  they  are  distributed  in  activation 
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energy  and  frequency  factor.  In  this  point  of  view,  which  we  have  taken  earlier 
[8,9],  in  a  narrow  range  of  observation,  the  relaxation  of  any  property  can  be 
characterized  by  a  simple  exponential  given  by: 


b 


% 


0i(O  = 


Xj{t)  -  Xj(oo) 
Xi{o)  -  X,(oo) 


(3) 


where  rpi(t)  is  the  level  of  the  normalized  property  at  time  t  of  an  element  t 
and  Xi(t),  Xj(o),  and  X,(oo)  are  the  actual  values  of  the  property  at  time  t, 
at  the  beginning  of  relaxation  and  at  the  end  of  relaxation,  respectively,  when 
equilibrium  is  reached.  Since  the  relaxation  is  a  thermally  activated  process,  it 
is  expected  to  have  a  characteristic  relaxation  time  r,,  given  by: 


n 


t,0  exp 


(4) 


where  A  Hi  is  an  activation  enthalpy  over  the  key  barrier  holding  back  the  local 
relaxation  and  l/r,0  is  a  fundamental  frequency  factor.  It  is  then  stated  further 
that  if  the  observation  is  not  in  a  narrow  range,  but  spans  over  many  decades  in 
time,  or  alternatively  is  accomplished  over  a  span  in  temperature,  other,  both 
easier  and  more  difficult  relaxation  processes  with  shorter  or  longer  relaxation 
times  become  also  observable.  This  requires  a  superposition  of  many  processes 
that  are  continuously  distributed  in  relaxation  time  r,  according  to  a  charac- 
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teristic  and  structure  dependent  normalized  distribution  function  p(r),  to  result 


in 


0(0=/  p(r)exp(-t/r)dr, 

Jo 


(5) 


and  obeying  the  normalization  condition 


(6) 


Many  procedures  of  operational  inversion  have  been  advanced  [10-13]  for  the 
determination  of  such  distribution  functions  from  experimental  data  collected 
over  long  periods  of  time  or  over  different  temperatures.  On  the  other  hand, 
based  on  quite  successful  empirical  procedures  dating  back  to  Kohlrausch  [14], 
more  than  a  century  ago,  a  different  approach  has  been  suggested  to  the  ac¬ 
counting  for  the  non-exponential  form  of  the  long  term  relaxation  behavior  of 
amorphous  media.  In  this  approach,  recently  revived  by  Williams  and  Watts 
[15],  it  is  observed  that  the  relaxation  in  the  normalized  property  can  be  fitted 
quite  well  over  a  relatively  wide  range  of  time  to  a  modified  relaxation  function 
with  a  fractional  exponent  given  as: 


where  tw  is  a  new  relaxation  time  and  /?  is  usually  in  the  range  of  0.5.  While 
many  recent  investigators  have  adopted  this  form  of  presentation  on  the  basis  of 
its  simplicity,  and  have  asserted  that  it  reflects  a  more  complex  non-Arrhenian 
relaxation,  others  have  searched  for  a  common  denominator  between  the  physi¬ 
cally  more  appealing,  but  operationally  more  complex  form  given  above  in  Eqn. 
(5)  and  the  simpler  form  given  by  Eqn.  (7)  [16,17].  Thus,  the  stretched  ex¬ 
ponential  form  of  Kohlrausch  has  been  found  consistent  with  a  random  walk 
of  thermally  activated  diffusion  of  some  impotent  configurations  to  a  central 
site,  where  upon  arrival,  a  measurable  unit  relaxation  event  takes  place.  The 
fractional  exponent  then  arises  from  a  distribution  of  free  energy  barriers  in  the 
preparatory  diffusion  process  of  the  key  impotent  configurations  to  the  central 
site  where  they  trigger  the  unit  relaxation  process.  To  obtain  the  fractional 
power  /?  in  the  exponent  of  the  exponential,  it  has  been  necessary  to  assume 
special  distributions  of  barrier  energies  with  exponential  tails.  Operational  pro¬ 
cedures  have  been  developed  further  to  obtain  distribution  functions  of  barrier 
heights  in  the  diffusion  problems  from  experimental  data  directly  [18]. 

We  consider  these  two  different  points  of  view  as  alternative  operational 
approaches  to  represent  the  same  physical  phenomenon,  in  which  the  relaxation 
process  is  distributed  in  relaxation  time  or  in  activation  energy  and  frequency 
factor.  Whether  the  actual  unit  relaxation  process  is  a  sudden  transformation 
of  a  local  configuration  from  a  more  open  one  to  a  more  compact  one  in  one 
thermal  fluctuation,  or  is  composed  of  a  series  of  prior  impotent  diffusive  steps 


approaching  a  central  site  without  any  observable  effect,  finally  triggering  the 
observable  unit  relaxation  process  upon  arrival  at  this  site  is  difficult  to  decide. 
The  topological  details  of  the  simulation  described  in  some  detail  in  (II)  [5] 
indicate  that  in  reality,  the  processes  may  be  a  combination  of  both  pictures, 
where  the  central  unit  process  is  noted  to  be  the  dissolution  of  5  and  7  sided 
polygon  dipoles.  Noting  difficulties  in  relating  the  Williams- Watts  relaxation 
function  uniquely  to  distribution  functions  of  energies  in  preparatory  diffusive 
processes  [18],  or  performing  an  operational  inversion  to  calculate  p(r)  from 
equations  of  the  type  given  by  Eqn.  (5),  we  will  merely  adopt  the  Williams- Watts 
function  and  obtain  the  best  fit  of  the  results  of  our  simulation  to  it.  We  justify 
this  less  than  complete  analysis  by  noting  that  the  simulation  covers  only  a  very 
small  portion  of  real  time  extending  over  only  200  atomic  fluctuations,  making 
it  virtually  impossible  to  obtain  any  relaxation  times  that  are  comparable  with 
experimental  information. 

The  results  of  fitting  the  ensemble  averages  of  the  evolving  structure  pa¬ 
rameters  presented  in  Fig.  la-lg  to  the  Williams-Watts  function  are  given  in 
Figs.  3a-3g.  In  these  figures,  the  overall  ensemble  averages  are  shown  by  the 
fluctuating  form.  The  best  fit  to  the  Williams-Watts  function  is  given  by  the 
smooth  solid  curve  to  determine  the  exponent  0  and  time  constant  tw.  Finally, 
the  dotted  curves,  also  shown  on  many  of  the  plots,  represent  the  best  fit  to  the 
simple  relaxation  function  with  0  =  1  (Eqn.  3).  The  resulting  best  fit  constants 
for  the  Williams-Watts  relaxation  function  are  given  in  Table  II.  We  note  from 
here  that  the  relaxation  of  atomic  volume,  enthalpy,  and  site  distortion  are  all 
characterizable  with  nearly  the  same  function  having  a  relaxation  time  in  the 
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vicinity  of  7  fluctuations  (=  140  time  steps)  and  an  exponent  of  /?  =  0.5.  The  re¬ 
laxation  of  the  shear  modulus  obeys  nearly  the  same  form  of  the  equation  with, 
however,  a  relaxation  time  of  more  than  twice  the  length.  In  comparison,  the 
relaxation  of  atomic  site  pressures  and  maximum  shear  stress  are  much  slower, 
as  is  apparent  from  Fig.  1,  and  is  characterizable  by  a  very  much  longer  time 
constant  of  around  53  fluctuations  (=  1060  time  steps)  and  a  considerably  larger 
exponent  of  /?  =  0.85  coming  close  to  the  simple  relaxation  process  of  Eqn.  (3). 

Examination  of  Figs.  3a-3g  show  that  all  ensemble  properties  decrease  from 
quenched-in  initial  excess  values  to  lower  values  during  the  simulation  span  of  200 
fluctuations.  The  decreases  are  more  dramatic  and  substantial  for  volume  per 
atom,  site  distortion  parameter,  and  enthalpy  than  they  are  for  the  individual 
stresses.  The  ensemble  averages  of  atomic  site  shear  moduli  actually  increase, 
as  they  should.  The  decrease  in  the  atomic  site  pressure  is  3%,  while  that  in 
atomic  site  maximum  shear  stress  is  5%.  The  related  rise  in  the  atomic  site 
shear  modulus  is  5%,  while  the  rise  in  the  atomic  site  bulk  modulus  could  be 
determined  to  be  6%  from  the  information  given  in  Fig.  lg.  The  corresponding 
decreases  in  the  ensemble  average  enthalpy  is  14%,  that  for  the  volume  per  atom 
is  10%,  and  for  the  site  distortion  parameter  only  2%.  Considering  that  the 
ensemble  average  fractional  change  of  enthalpy  must  be  predominantly  coming 
from  the  change  in  internal  strain  energy  stored  partly  in  shear  strain  energy 
and  partly  in  bulk  compression  strain  energy,  we  estimate  the  fractional  enthalpy 
change  from 


2 

9* 


(8) 


La*  A. 


A/i  a  /i  f  P2  T2max\ 

T~Aln[2K  +  ~^r) n’ 


where  K  is  the  bulk  modulus  and  /z  the  shear  modulus.  Evaluations  of  Eqn. 
(8)  gives  for  Ah/h  exactly  14  %  from  the  changes  in  pressure,  maximum  shear 
stress,  bulk  modulus,  shear  modulus,  and  atomic  volume.  We  consider  this  as  a 
good  check  on  internal  consistency  in  the  simulation. 

We  note  further,  in  passing,  that  the  ratio  of  the  ensemble  average  max¬ 
imum  shear  stress  to  shear  modulus  is  0.21,  while  the  corresponding  ratio  of 
average  pressure  to  bulk  modulus  is  0.085.  These  are  both  somewhat  higher 
than  the  ratios  reported  by  Egami  and  Vitek  [19]  for  a  three-dimensional  amor¬ 
phous  medium.  We  attribute  the  difference  to  the  two-dimensional  nature  of 
our  material. 

Since  the  relaxation  of  most  of  the  properties,  except  the  atomic  site  stresses, 
can  be  characterized  by  a  Williams- Watts  function  with  a  fractional  exponent 
of  0.5  and  a  characteristic  time  constant  rw  =  7  atomic  fluctuations  (=  240  time 
steps),  it  is  possible  to  give  a  distribution  function  p(r)  of  relaxation  times  of 
simple  Arrhenian  processes  on  the  basis  of  an  analytical  solution  of  the  integral 
equation  obtained  by  equating  Eqns.  (5)  and  (7).  This  solution,  obtained  by 
Lindsay  and  Patterson  [18],  gives: 


pW  =  ib^exp(^) 


1 
m 

m!1 


For  comparison  with  experimental  results  obtained  by  Deng  and  Argon  [8], 
we  convert  this  distribution  into  one  of  activation  free  energies  AG*  by  noting 
that 

p  (AG*)  dAG*  =  p(r)dT  =  p  d  ,  (10) 

giving 

=  <n> 

Upon  introduction  of  the  relations 

(£)  =  © exp  (if)  ’  r<J  =  1,2l"'a-  (12“-6) 

into  Eqns.  (9)-(ll)  it  is  found  that: 
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where  i/c  is  a  fundamental  atomic  cluster  frequency  related  to  the  clusters  which 
undergo  structural  relaxations.  Noting  that  for  the  simulation  kT/E0  =  0.1,  and 
taking  Ug  as  the  atomic  frequency,  tq/tw  =  l/l 4t  for  tu  =  7tq,  as  obtained  from 
Table  II,  it  is  possible  to  evaluate  the  distribution  function  of  Eqn.  (13)  as  a 
function  of  A G*  jEa.  The  result  of  this  is  given  in  Fig.  4.  While  the  shape  of 
this  figure  compares  favorably  with  that  obtained  by  Deng  and  Argon  [9]  from 
internal  friction  measurements  in  a  Cu^Zr^i  alloy  (their  Fig.  8),  the  energy 
scale  of  the  latter  is  a  factor  of  26.8  larger  than  the  simulation  results.  This  is 
not  surprising,  since  the  frequencies  in  the  internal  friction  results  were  in  the 
range  of  0.25  HM,  while  those  in  the  simulation  were  no  smaller  than  lO11#*. 
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3.3  Inhomogeneities  and  Clustering  of  Atomic  Motions 


Figures  5a  and  5b  show  the  total  atom  motions  in  two  separate  time  periods 
of  10  fluctuations  each  (5a  from  70-80;  5b  from  100-110).  The  displacements,  at 
a  magnification  of  10,  are  shown  as  bars  emanating  from  the  initial  positions  of 
the  atoms  at  the  beginning  of  the  time  period.  The  motions  shown  in  these  two 
frames  are  typical  of  many  such  frames  that  have  been  examined.  In  Fig.  5a, 
the  motions  are  relatively  random  over  most  of  the  field,  except  in  the  central 
upper  portion,  where  group-like  displacements  toward  the  left  are  evident.  In 
Fig.  5b,  on  the  other  hand,  very  large  group-like  motions  of  atoms  toward  the 
left  are  visible  in  the  entire  right  half  of  the  simulation  cell.  Such  large  to  and  fro 
motions  of  blocks  of  atoms  often  reversed  in  the  following  10  fluctuations.  They 
are  attributed  to  the  resonances  resulting  from  the  boundary  mass  that  show 
up  clearly  in  the  pulsations  of  the  volume  shown  in  Fig.  la  and  were  therefore, 
largely  ignored. 

More  importance  was  attached  to  the  clustering  of  the  atomic  structure  pa¬ 
rameters,  which  were  considered  to  be  more  genuine  and  were  followed  in  some 
detail.  Figures  6a-6c  show  the  clustering  of  enthalpy,  volume  per  atom,  and  site 
distortion  parameter  averaged  over  the  entire  simulation  span  of  200  fluctua¬ 
tions.  The  length  of  the  solid  horizontal  bars  reflect  the  relative  magnitudes  of 
the  respective  properties.  All  the  atomic  site  properties  that  have  been  encir¬ 
cled,  are  larger  than  the  ensemble  average,  except  those  in  circles  3  and  5,  which 
are  smaller  than  the  ensemble  average.  The  bars  emanate  from  the  centers  of 
atoms.  The  atoms  bearing  numbers  are  thus  either  in  the  upper  or  lower  10 


percentile  of  the  population  with  largest  and  smallest  properti  Clustering  in 
atomic  level  pressure,  maximum  shear  stress,  xy  shear  stress,  and  shear  modu¬ 
lus,  or  bulk  modulus  was  less  prominent.  For  these,  the  excess  properties  which 
were  also  large  were,  however  more  uniformly  mixed.  These  clusters  in  enthalpy, 
volume  per  atom,  and  site  distortion  indicate  that  excess  properties  in  the  field 
have  often  long  life  times  and  persist  through  even  relatively  long  simulation 
spans.  If  they  are  to  represent  actual  structural  relaxations  in  real  time,  they 
must,  of  course,  have  life  times  that  must  be  larger  by  orders  of  magnitude  than 
the  period  of  200  fluctuations  considered  here.  The  evolution  of  properties  in 
clusters  and  their  random  changes  in  size  and  shape  was  followed  in  more  de¬ 
tail  over  smaller  increments  of  the  simulation  span,  but  maps  of  these  changes 
will  not  be  presented  here.  From  examination  of  these  changes,  however,  the 
following  important  observations  have  been  made: 

a.  Prominent  clustering  occurs  only  in  excess  enthalpy,  volume  per  atom  (free 
volume),  and  site  distortion.  These  are  also  the  properties  that,  as  shown 
in  Figs.  1  and  2  show  the  more  prominent  relaxation  effects. 

b.  Clustering  of  atomic  level  stress  is  slight.  The  largest  clusters  of  atoms 
with  excess  pressure  or  excess  maximum  shear  stress  were  only  2-3  atoms 
across.  Generally,  it  was  found  that  atoms  under  large  stress  of  all  types, 
i.e.,  pressure,  or  maximum  shear  stress  were  the  Zr  atoms,  while  those 
under  the  lowest  pressure  or  shear  stress  were  the  Cu  atoms.  Since  these 
were  initially  randomly  spaced,  and  since  the  simulation  time  spans  were 
too  short  to  produce  much  chemical  short  range  ordering,  the  lack  of  clus¬ 
tering  is  not  surprising.  As  already  noted,  the  clustering  correlates  with 
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more  prominent  relaxation  strength.  Somewhat  similar  observations  were 
made  also  by  Maede  and  Takeuchi  [20]  in  their  study  of  radial  correlation 
functions  of  the  atomistic  structure  parameters,  but  they  did  not  give  a 
cause. 

c.  Comparison  of  the  clustering  in  Figs.  6a-6c  shows  that  the  clusters  with 
high  excess  enthalpy  tend  to  overlap  with  clusters  having  high  excesses  of 
free  volume,  and  those  showing  excesses  of  site  distortior. .  Parenthetically, 
it  must  be  observed  that  clustering  of  excess  enthalpy  requires  clustering 
of  a  combination  of  dilatational  and  shear  strain  energy,  which,  however, 
does  not  require  clustering  of  either  pressure  or  maximum  shear  stress. 

d.  Relief  of  site  distortion  in  clusters  of  this  property  resulted  in  clustering 
of  shear  strains  of  random  sign  shown  in  Fig.  6d,  which  showed  spatial 
correlation  with  the  clusters  of  site  distortion.  A  similar  correlation  was 
found  between  the  relief  of  excess  free  volume  and  volumetric  compaction. 

e.  Some  spatial  correlation  was  found  between  local  shear  strain  production 
during  relaxation  and  increments  in  positive  dilatation.  Such  dilatancy 
is  known  to  be  a  prominent  effect  during  imposed  shear  strains.  The 
importance  of  this  effect  is  more  clearly  demonstrated  in  the  comparison 
study  (IV)  of  simulation  of  plastic  shear  deformation  [16]. 

f.  It  was  found  generally  that  regions  with  excesses  of  free  volume  coincided 
with  regions  with  reduced  bulk  modulus,  or  shear  modulus,  or  both,  while 
the  opposite  was  true  for  regions  of  deficiency  of  free  volume.  This  is 
considered  to  be  an  important  source  of  distribution  in  frequency  factor 
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found  by  Deng  and  Argon  [9]  in  internal  friction  experiments  in  CuS9Zr4l 
glasses. 

g.  Comparison  of  the  clustering  of  site  distortion  with  distributions  of  Voronoi 
polygons  with  sides  other  than  six  reported  in  the  companion  study  (II) 
(5],  showed  a  necessary  correlation  between  large  distortions  and  non- 
hexagonal  polygons. 

In  summary  of  the  above  observations,  we  note  that  three  types  of  clusters 
were  observed:  an  excess  volume  cluster  (free  volume),  called  an  n  cluster  by 
Egami  and  Vitek  [19],  has  a  lower  local  pressure  and  a  lower  local  modulus;  a 
cluster  of  deficiency  of  volume,  called  a  p  cluster  by  Egami  and  Vitek,  which 
possesses  a  higher  local  pressure  and  a  higher  local  modulus;  and  finally,  a  site 
distortion  cluster,  which  possesses  a  higher  level  of  maximum  shear  stress  and 
also  an  excess  of  volume.  Only  the  first  and  third  types  of  cluster  undergo 
important  evolution  during  structural  relaxation,  while  the  first  contracts  the 
third  undergoes  a  combined  shear  and  dilatation.  The  changes  in  the  p  type 
cluster  are  far  less  consequential  or  apparent. 

The  evolution  of  properties  in  5  clusters  were  studied  in  more  detail  as  a 
function  of  time  during  the  simulation.  These,  marked  1-5  in  Figs.  6a  and  6c, 
each  contain  roughly  equal  numbers  of  Zr  and  Cu  atoms.  Cluster  1  possessed 
a  substantial  excess  of  volume  per  atom  in  the  initial  state.  Clusters  2  and  4 
possessed  a  large  excess  in  site  distortion,  while  clusters  3  and  5  had  quite  low 
levels  of  excess  volume  per  atom.  The  changes  in  volume  per  atom,  enthalpy, 
site  distortion,  atomic  site  pressure,  maximum  shear  stress,  atomic  site  shear 
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modulus,  are  shown  as  a  function  of  time  of  simulation  in  Figs.  7a-7f.  We  note 
from  examination  of  these  traces  that  cluster  1,  which  has  the  most  extreme 
excesses  in  volume  shows  the  most  dramatic  relaxations  in  volume,  in  enthalpy, 
in  site  distortion,  in  rise  in  pressure,  and  in  rise  in  shear  modulus,  mostly  by 
undergoing  a  decrease  in  dilatation.  Clusters  3  and  5  with  the  minimum  excesses 
in  properties,  show  no  important  changes,  but  exhibit  only  flutter.  Clusters  2 
and  4,  which  had  large  initial  site  distortions,  show  larger  variations  in  enthalpy, 
site  distortion,  pressure,  and  shear  stress,  but  no  systematic  changes. 

IV.  DISCUSSION 

Isothermal  structural  relaxations  in  real  amorphous  media  at  about  half  the 
glass  transition  temperature  would  be  comparatively  very  slow  processes.  In 
metallic  glasses  such  as  the  CuTZf\-t  glasses  that  have  been  simulated  here, 
structural  relaxations  have  characteristic  time  constants  that  are  closely  compa¬ 
rable  to  those  found  in  anelastic  relaxations  |9]  and  have  activation  energies  that 
span  between  1-2  eV  and  frequency  factors  that  are  in  the  range  of  1010  —  1012s-1. 
This  gives  ranges  in  relaxation  time  from  10s  —  10lss.  Processes  with  such  very 
long  relaxation  time  can,  of  course,  not  be  simulated  by  MD.  Therefore,  in  the 
simulation  reported  here,  rather  artificial  initial  conditions  were  taken,  in  which 
an  alloy  at  its  melting  point  was  quenched  at  infinite  rate  to  a  low  temperature, 
where  it  had  an  extremely  large  excess  of  properties  that  could  be  relaxed.  The 
time  span  in  the  relaxation  that  could  be  reasonably  simulated  consisted  of  only 
200  thermal  fluctuations  at  a  temperature  equal  to  0.556  of  Tg,  i.e.,  over  a  period 
differing  from  the  required  real  times  by  10-11  decades.  Consequently,  the  results 
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presented  here  cannot  be  compared  directly  with  those  of  metallic  glasses  under 
any  realistic  conditions.  Therefore,  these  results  must  be  examined  primarily  for 
the  qualitative  information  that  they  convey,  for  the  chronological  ordering  of 
processes,  and  the  time  independent  kinematics  of  structural  relaxation.  This  is 
also  apparent  from  the  distribution  function  of  activation  energies  for  relaxation 
of  excess  enthalpy,  volume  per  atom,  and  distortion  that  were  calculated  by 
matching  a  Williams- Watts  relaxation  function  to  the  results  of  the  simulation. 
The  distribution  function  has  a  peak  at  A G* / E0  =  0.45.  Since  the  scaling  factor 
E0  =  0.15eV.  for  this  simulation  [4],  the  observed  characteristic  activation  energy 
is  0.068  eV.  and  not  1.5  eV.,  as  is  the  case  for  the  real  CuxZri_x  alloy.  Clearly, 
the  only  processes  that  can  be  relaxed  over  such  short  periods  of  time  are  those 
that  appear  in  the  simulation.  Indeed,  if  this  activation  energy  and  its  associated 
relaxation  period  of  7  fluctuations  are  taken  as  accurate,  and  a  formal  upward 
scaling  of  relaxation  times  to  r2  =  4  seconds  (a  0.25  Hz  internal  friction  experi¬ 
ment)  at  T  =  0.556 Tg  of  a  Cu^ZrAi  glass  having  a  Tt  of  741K  is  done,  according 
to  direct  scaling  relation  of  A G\  =  A G\  +  kTlnfa/T])  an  activation  energy  of 
AGJ  =  1.05  eV  is  obtained,  where  AGJ  =  0.068eF,  T\  =  7x5.4xl0_13sec. 
and  T  =  412K  have  been  used  from  the  simulation.  This  is  almost  exactly  the 
lowest  energy  in  the  spectrum  of  energies  measured  in  internal  friction  experi¬ 
ments  by  (at  0.25  Hz)  by  Deng  and  Argon  [9],  that  can  produce  relaxation  at 
such  a  low  homologous  temperature.  Thus,  when  properly  scaled  up,  the  results 
of  the  simulation  are  consistent  also  with  experimental  results.  This  however, 
implies  a  similarity  I  a  the  shapes  of  the  energy  spectra  at  the  atomic  fluctuation 
scale  with  energies  in  the  macroscopic  scale  for  which  no  clear  justification  can 
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be  furnished. 


Considering  the  results  of  the  simulation  as  a  whole,  we  have  shown  that 
structural  relaxations  are  complex  and  cooperative  processes.  The  excess  struc¬ 
tural  properties  in  the  form  of  excess  enthalpy,  volume  per  atom  (free  volume), 
site  distortion,  atomic  level  stresses,  and  moduli  are  not  uniformly  spread  out 
in  the  structure.  The  excess  enthalpy,  free  volume,  and  site  distortion  show  the 
greatest  degree  of  clustering,  while  the  atomic  level  stresses  cluster  less,  since 
they  are  associated  more  directly  with  the  Zr  and  Cu  atoms,  which  were  ran¬ 
domly  distributed  in  the  initial  state  of  the  structure.  At  the  temperature  of 
simulation  (0.556T,),  the  relaxation  is  to  a  large  extent  a  local  and  immobile 
one,  i.e.,  the  clusters  with  large  excesses  in  properties  show  a  decay  in  these 
by  local  atom  motions.  There  were  no  important  displacements  of  regions  with 
large  free  volume  toward  regions  with  deficiencies  in  free  volume,  i.e.,  so-called 
n  type  clusters  were  not  observed  in  any  meaningful  way  to  migrate  to  p  type 
clusters,  as  described  by  Egami  and  Vitek  [19]  in  their  simulations.  While  such 
“interdiffussion”  of  defect  clusters  is  very  appealing,  and  may  probably  occur 
in  more  realistic  situations  in  real  time,  they  must  require  much  longer  periods 
of  time,  in  excess  of  the  simulation  period  studied  here.  We  recall  from  the 
associated  study  (II  [5]),  that  the  topological  features  of  the  relaxation  involve 
a  reduction  in  the  volume  fraction  of  atom  sites  with  5  and  7  near  neighbors, 
i.e.,  reduction  of  liquid-like  boundary  material  by  gradual  dismemberment  and 
dissipation  of  such  regions,  as  they  are  replaced  by  a  growth  of  quasi-ordered 
material  consisting  of  atoms  with  6  near  neighbors.  Our  last  associated  study 
(IV)  [6]  demonstrated  that  it  is  this  liquid-like  boundary  material,  which  seeds 


24 


jwaawjwj  'A  iw.v .vw-v.  >■'  r.  >■  re  j.  >.-j  ij  ? .  y  vv  t..  ^  -.  <.  -v*v<- 


shear  transformations  that  are  the  principal  modes  of  strain  production  in  aged 
amorphous  material  before  the  ordered  domains  become  large  enough  to  permit 
plasticity  by  dislocation  glide. 
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TABLE  II.  Best  Fit  Values  of  Williams- Watts  Distribution  Parame¬ 
ters  Determined  from  Relaxation  Simulation 


Structure  Parameter 

fi 

W 

*0 

*oo 

Atomic  Volume 

0.55 

130 

1.668 

1.504 

Enthalpy 

0.5 

140 

-1.344 

1.565 

Distortion  parameter 

0.55 

150 

1.082 

1.060 

Atomic  Pressure 

0.85 

1060 

1.486 

1.454 

At.  Max.  Shear  Stress 

0.85 

1020 

0.968 

0.916 

Shear  Modulus 

0.5 

320 

4.232 

4.396 

(*)  in  time  steps  (=0.05  fluctuations) 
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APPENDIX  I 


Boundary  Masses 

To  maintain  the  external  pressure  constant  in  MD  simulations,  Andersen  [7] 
introduced  another  dynamic  coordinate  into  the  basic  Lagrangian  representing 
the  equations  of  motion  of  the  atoms  in  the  simulation.  Th  new  coordinate 
Q  representing  the  volume  of  the  simulation  cell  then  gives  the  modified  La¬ 
grangian: 


£(*,*,<?,<?)  =  -\'h<t>{Q1,zPa)  +  \MQ2-(*Q  (^-i) 


i=i 


where  all  quantities  have  been  defined  earlier  in  Section  2.2  of  the  text. 

To  minimize  artifacts  of  unwanted  resonances  of  the  fictitions  boundary 
masses  and  the  atoms  in  the  simulation  cell,  it  is  necessary  to  choose  the  bound¬ 
ary  mass  carefully,  not  to  affect  the  kinetics  of  any  process  that  is  being  simu¬ 
lated.  The  recommended  procedure  is  to  chose  M  in  such  a  way  that  the  period 
of  the  fluctuations  of  Q  is  approximately  equal  to  Q1^3  divided  by  the  speed 
of  sound  in  the  simulation  cell.  As  discussed  in  Section  2.2,  this  criterion  gave 
a  boundary  mass  of  4  atom  masses.  Nevertheless,  to  explore  the  effect  of  the 
boundary  mass  on  the  results  of  the  simulation,  other  values  one  hundredth  of 
the  chosen  value  and  hundred  times  the  chosen  value  were  also  used  to  determine 
the  effect  of  the  choice. 
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Figures  A-la  -  A-lc  show  the  results  of  a  relaxation  simulation  of  the  volume, 
enthalpy,  and  internal  pressure  for  the  three  different  boundary  masses.  Clearly, 
the  responses  to  boundary  masses  of  0.04  and  400  give  insufficient  damping 
and  erratic  behavior  or  supercritical  damping  respectively  during  the  period 
of  simulation  lasting  1000  time  steps  (50  fluctuations).  The  results  for  the 
boundary  mass  of  4  show  rather  regular  oscillations  in  properties  that  were 
presented  already  above.  These  could  be  readily  isolated  and  discounted  when 
necessary. 
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FIGURE  CAPTIONS 


Fig.  1 


Fig.  2 


Fig.  3 

Fig.  4 


Fig.  5 


Fig.  6 


Time  dependent  changes  in  ensemble  averages  over  200  fluctuations  (= 
4000  time  steps):  (a)  volume  per  atom,  (b)  enthalpy,  (c)  site  distortion 
parameter,  (d)  atomic  level  pressure,  (e)  atomic  level  maximum  shear 
stress,  (f)  xy  shear  stress,  (g)  atomic  site  bulk  modulus.  Ensemble 
averages  are  given  both  for  the  entire  population  as  well  as  for  the  Zr 
and  Cu  atoms  separately. 

Time  dependent  changes  in  root  mean  square  deviations  from  the 
ensemble  mean  values  over  200  fluctuations:  (a)  volume  per  atom,  (b) 
enthalpy,  (c)  site  distortion  parameter,  (d)  atomic  level  pressure,  (e) 
atomic  level  maximum  shear  stress,  (f)  xy  shear  stress,  (g)  atomic  site 
bulk  modulus. 

Best  Williams- Watts  function  fits  to  the  changes  in  ensemble  averages 
for  the  same  seven  structure  parameters  given  in  Figs.  1  and  2. 

The  free  energy  barrier  distribution  function  for  the  structural  relax¬ 
ations,  obtained  from  the  Williams- Watts  distribution  function  by  an 
operational  inversion. 

Two  different  states  of  incremental  atom  displacements  during  the 
simulation:  (a)  for  the  period  between  70-80  fluctuations,  (b)  for  the 
period  between  100-110  fluctuations. 

Field  of  excess  properties  averaged  over  the  entire  simulation  period: 
(a)  volume  per  atom,  (b)  enthalpy,  (c)  site  distortion  parameter,  (d) 
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relaxation  shear  strains. 

Fig.  7  Time  dependent  changes  of  the  structure  parameters  in  five  selected 
clusters:  (a)  volume  per  atom,  (b)  enthalpy,  (c)  site  distortion,  (d) 
atomic  site  pressure,  (e)  atomic  site  maximum  shear  stress,  (f)  atomic 
site  shear  modulus. 

Fig.  A1  Effect  of  choice  of  boundary  mass  on  the  relaxation  of  excess  properties: 
(a)  volume  per  atom,  (b)  enthalpy,  (c)  atomic  site  pressure. 
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Time  dependent  changes  in  ensemble  averages  over  200  fluctuations  (= 
4000  time  steps):  (a)  volume  per  atom,  (b)  enthalpy,  (c)  site  distortion 
parameter,  (d)  atomic  level  pressure,  (e)  atomic  level  maximum  shear 
stress,  (f)  xy  shear  stress,  (g)  atomic  site  bulk  modulus.  Ensemble 
averages  are  given  both  for  the  entire  population  as  well  as  for  the  Zr 
and  Cu  atoms  separately. 
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Fig.  2 


Time  dependent  changes  in  root  mean  square  deviations  from  the 
ensemble  mean  values  over  200  fluctuations:  (a)  volume  per  atom,  (b) 
enthalpy,  (c)  site  distortion  parameter,  (d)  atomic  level  pressure,  (e) 
atomic  level  maximum  shear  stress,  (f)  xy  shear  stress,  (g)  atomic  site 
bulk  modulus. 
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Fig.  5  Two  different  states  of  incremental  atom  displacements  during  the 
simulation:  (a)  for  the  period  between  70-80  fluctuations,  (b)  for  the 
period  between  100-110  fluctuations. 
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Time  dependent  changes  of  the  structure  parameters  in  five  selected 
clusters:  (a)  volume  per  atom,  (b)  enthalpy,  (c)  site  distortion,  (d) 
atomic  site  pressure,  (e)  atomic  site  maximum  shear  stress,  (f)  atomic 
site  shear  modulus. 
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Fig.  A1  Effect  of  choice  of  boundary  mass  on  the  relaxation  of  excess  properties: 
(a)  volume  per  atom,  (b)  enthalpy,  (c)  atomic  site  pressure. 
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